Variational study of the extended Hubbard-Holstein model on clusters of variable site 
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We study the complete extended Hubbard-Holstein Hamiltonian on a four-site chain with equally 
spaced sites, with spacing-dependent electronic interaction parameters evaluated in terms of Wannier 
functions built from Gaussian atomic orbitals. By successive application of generalized displacement 
and squeezing transformations, an effective polaronic Hamiltonian is obtained, containing a purely 
phonon-induced, long-range intersite charge interaction. The phase diagrams for A'^ — 1,2,3,4 
electrons have been determined by variational minimization of the sum of the electronic and the 
phononic energy, paying special attention to the effects of the above mentioned phonon-induced 
long-range interaction. To characterize the physics of each ground state, we have evaluated how 
variations of the site spacing affect the behavior of several correlation functions of experimental 
interest, with the aim of representing lattice deformation effects, either spontaneous or induced by 
external pressure. 



I. INTRODUCTION. 

The model describing correlated electrons in a narrow 
band interacting locally with lattice deformations, the 
so-called Holstein-Hubbard model, has been studied in- 
tensively in the last years both on infinite lattices ||] 
and on small clusters jl^- ||l8|. The classes of materials 
to which the model is assumed to apply, such as, for in- 
stance, the high-temperature superconducting cuprates 
or the manganites showing colossal magnctorcsistance, 
are characterized by rather strong electron-phonon cou- 
pling, so that a perturbative treatment is in general un- 
reliable. In that case, one can apply to the interacting 
electron-phonon Hamiltonian a unitary "displacement" 
transformation [ p9| generated by a suitably chosen anti- 
hermitian operator D. From the displaced Hamiltonian 
H = H a non-perturbative effective Hamiltonian 
for phonon-dressed electrons (polarons) can be obtained 
by averaging H over an appropriate phonon wavefunc- 
tion. In recent years, the use of variationally optimized 
"squeezed" phonon states of the form e^^\Qph), with 
\Oph) being the phonon vacuum, has proved to be advan- 
tageous ^,^,|ll|,0 over methods using simple harmonic 
oscillator states. However, as the operators D and S on 
a lattice with L sites each contain L undetermined c- 
numbers, one for each phonon wavevector q, such proce- 
dure requires to set all such 2L displacement and squeez- 
ing parameters to their optimal values. A variational de- 
termination on an infinite lattice is therefore impossible, 
and one has to adopt perturbative techniques |2l|,^ . For 
this reason, the most popular approach has been to ne- 
glect the (/-dependence of the displacement parameters, 
by assuming 6q = 1 for any q. This is the standard Lang- 
Firsov approximation (LFA) which, however, has at 
least two limitations. On one side, it neglects any inter- 



site phonon correlation, preventing the study of how the 
charge on a given site can induce deformations also on 
different sites, possibly giving rise to extended polarons 
[Q. On the other side, one can show |^ that, for non- 
dispersive phonons, the LFA leads to the vanishing of a 
purely phonon-induced intersite charge interaction (the 
"residual interaction" first introduced in Ref. ||2^), in 
principle present in H. One concludes that the physics 
of the system is not faithfully represented in LFA. If, in- 
stead of infinite lattices, one considers small clusters [ p^ , 
the numerical complexity is reduced, and one can rigor- 
ously take into account the g-dependence of all the varia- 
tional parameters. In this way the limitations of LFA can 
be avoided, and, in particular, the inter-site correlations 
between charges and deformations, as well as the effects 
of the purely phonon-induced interaction terms, can be 
satisfactorily taken into account. 

This is actually one of the points on which the present 
work focuses. Indeed, we shall first derive the effective 
polaronic Hamiltonian from a bare Hamiltonian which 
includes in the electronic part all the one- and two-body 
interactions compatible with a non-degenerate band of 
electrons, and has an electron-phonon coupling of the 
Holstein type. As far as the set of electronic interactions 
is concerned, the present work extends previous stud- 
ies, where the only two-body electronic coupling terms 
included in the bare Hamiltonian were the on-site and 
inter-site Coulomb charge interactions. Then we shall 
solve exactly the effective model on a four-site chain with 
constant intersite distance a and closed boundary condi- 
tions. Considering that in this case the wavevector q 
takes four values, given by — 7r/2a, 0, 7r/2a and n/a, the 
independent displacement and squeezing parameters are 



six in total, as S, 



Tr/2a 



S-7r/2a and a^/2a 



-7r/2a 



They 



will be simultaneously and independently determined by 



1 



exact diagonalization of the effective Hamiltonian, fol- 
lowed by variational optimization of the ground state en- 
ergy, for given, but arbitrary, site spacing and given shape 
of the Wannier functions, as discussed below. This will 
allow us, on the one hand, to discuss the interplay be- 
tween charges and phonons both on the same site and on 
different sites, and, on the other hand, to properly take 
into account the effect of the phonon-induced inter-site 
charge interaction. 

Another novel aspect of the present work is that, as 
we model the dependence on the site spacing of all the 
electronic interaction parameters, we can study the phase 
diagram and the correlation functions at any given fill- 
ing as functions of the intersite distance. Our results 
can therefore give qualitative indications relevant both 
to measurements under mechanical or chemical pressure, 
and to materials showing "stripes" ||2^ or other forms 
of spontaneously coupled charge and structural inhomo- 
geneities , where regions with different lattice param- 
eters, and different electronic character, have been re- 
ported to coexist. 

The paper is organized as follows. In Section 2, we de- 
fine the model and derive the effective polaronic Hamil- 
tonian. The results of the numerical analysis on a four- 
site closed chain for different fillings N = 1,2,3,4 are 
discussed in Section 3. Section 4 is devoted to the con- 
clusions. 

II. THE MODEL AND THE EFFECTIVE 
HAMILTONIAN. 

Our model Hamiltonian H = H^i + Hei-ph + Hph con- 
sists of a generalized Hubbard Hamiltonian for correlated 
electrons coupled to phonons of frequency Qq. The bare 
electronic Hamiltonian reads 

Hei = ej^ja - ^[tji - X^i{n^ _„ + n,__^)]c]^c,^ 

3 1 

where cj^ (cjo-) creates (annihilates) an electron with spin 
a on site j, and Sj and rij are the total spin and charge on 
site J, respectively. Following the procedure reported in 
Ref. the bare electronic interactions e,t,U,V, J, X , 
are parametrized by means of Wannier function built 
from Gaussian functions (f>j{r) — (2r^/7r)'^/^ exp[— r^(r — 
Rj)'^] modelling the atomic orbitals on the j-th site in the 
position Rj. The Wannier functions, and by consequence 
the interactions, depend on the intersite distance a and 
the width F"^ of the orbitals. The explicit expressions 
of e,t,U,V, J, P, X that we use in this paper have been 



evaluated in the case of a dimer and are given in Ref. [g^ . 
This allows us to estimate all the electronic parameters 
in an ab-initio type of scheme. 

The interaction Hei-ph between the electrons and the 
phonons is assumed of the Holstein type Hei-ph — 
G^j^UjUjc;, where Uj is the lattice deformation at 
site j. By introducing the characteristic length Cq = 
^Jfi/2MVlq and the renormalized electron-phonon cou- 
pling Qq = GCq, in second quantization Hgi-ph takes the 
form 

Hei-ph ^Y^i (^-9 + ^9) "9 
9 

where Uq = L~^/^ ^ Uj ^ exp[iqRj], L being the num- 
ber of lattice sites. Finally, the free phonon term 

Hp,, = Y^nqiblbq + 1/2) (3) 
9 

completes the model Hamiltonian. 

The general procedure to obtain the effective polaronic 
Hamiltonian and to identify its ground state consists in 
the following steps : 

(i) a unitary transformation generated by a suitably 
chosen anti-hermitian operator D is used to map the 
model Hamiltonian H onto the "displaced" Hamiltonian 
H = e^He-° |l9). We choose 

^ = E'^.^(^-9-M"9 (4) 
9 

where the parameters 5^, defined on the L wavevectors 
of the first Brillouin zone, measure the degree of "dis- 
placement" of the corresponding phonon mode. Their 
real-space expressions 5ji^ associated with the coupling 
between charge on site j and phonons on site I , are given 
by 5ji = (1/L) (5g exp[ig(i?j — Ri)]. It is evident that 
assuming 5q ~ 1 for any q (as it is done in the LFA) im- 
plies 5ji = if Rj ^ Ri, forcing the effective interactions 
to be strictly local; 

(ii) by averaging H over the squeezed phonon state 
l^p/i) 1^1^ the phononic degrees of freedom are elimi- 
nated, yielding the effective polaronic Hamiltonian. The 
squeezed state is defined as \'^ph) = e~^\Oph) 0,^,^, 
where \Oph) is the harmonic oscillator vacuum state, such 
that bq\Oph) — 0, and the squeezing generator S, ex- 
pressed in terms of the L variational parameters aq (or, 
in real space, aji = (l/L) J2q <^q e:xjp[iq{Rj — Ri)]) reads: 

^ - E (^9^-9 - ^9^-9) = E - ^^^0 (5) 
9 j,i 

This procedure leads to an effective Hamiltonian H* de- 
scribing phonon-dressed electrons (which we shall gener- 
ally designate as polarons) having the form 
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sinh^ (ofp) + 



jla 



with the renormalized interactions 

J2 



(6) 



(7) 



n are assumed, then this phonon-induced long-range 
part of the Coulomb interaction remains non-vanishing 
only if the ^-dependence of the displacement parame- 
ters is explicitly accounted for. This interaction was 
first discussed within LFA, for dispersive phonons, in 
ref. p6[ , but with no attempt of estimating its value. 
A more recent variational study employing a non- 
factorized fermion-boson wavefunction in the context of 
a simple Hubbard electronic Hamiltonian (i.e., such that 
Xij ~ Vij = Jij = Pij = 0), also found an analogous 
effective interaction. 



III. RESULTS 



(8) 



n = ^.-^ - i E ^ (2 - ^''^'''-'''^ (9) 



t*i = Ttji X*i = r Xji P*i = t'^Pji Jji = Jji 



(10) 



exp 




cos (qa)] e 



-2q„ 



(11) 



(iii) the eigenstates and eigenvalues of the polaronic 
Hamiltonian are determined by using an exact diagonal- 
ization method proposed in Ref. [Q; 

(iv) the optimal values of the displacement and squeez- 
ing parameters are determined variationally by minimiz- 
ing the total energy of the electron-phonon system. 

A detailed derivation of Eqs. (0-|ll|) is given in Ref. . 
We see that all the parameters, but J, are renormal- 
ized to some extent by the phonons. The one- and two- 
particle hopping amplitudes (respectively t*i , X*i and 
P*i) depend on the hopping reduction factor r, whose 
value, once the optimal Sq and Uq are given, follows from 
Eq. (ITll) . The charge mobility is thus controlled simulta- 
neously, and competitively, by the displacement and the 
squeezing parameters. 

In the following, all the bare intersite coupling con- 
stants will be assumed non-vanishing only between 
nearest-neighbor sites. We stress that, even under this 
assumption, the effective intersite density-density inter- 
action V*i (Eq.(9)) acquires a phonon-dependent long- 
range contribution generally extending beyond nearest- 
neighbor sites. If dispersionless phonons of frequency 



The general procedure illustrated in Section 2 has been 
applied to the case of a four-site chain with dispersionless 
phonons of frequency Tifl = 0.1 eV, a value appropriate 
to cuprates and manganites. The eigenvalues and the 
eigenvectors of the effective Hamiltonian have been ex- 
actly calculated using the method outlined in Ref. Q. 
The knowledge of all the eigenvalues for given a and T al- 
lows to identify the ground state by the independent and 
simultaneous optimization of the wavevector-dependent 
variational parameters defining, for each phonon mode 
q, the degree of displacement (Sq) and squeezing (ctg). 
Once the optimal values of the parameters are identified, 
the full eigenvalue spectrum can be obtained. The value 
of r might also have been obtained variationally, as in 
Ref. ||l^,|l^ , but due to the increased complexity of the 
numerical analysis we have investigated the phase dia- 
gram only for three representative values of F, equal to 
0.8, 1, 1.2 A"^ We find that T'^ roughly scales as a 
in determining the main features of the phase diagram. 
The differences in the phase diagram for different T val- 
ues at given N are only quantitative, and correspond to 
slight variations of the border lines between the different 
phases. In particular, we have found that for increasing 
N the effects induced by variations of F tend to become 
less and less significative. Therefore, we will limit our- 
selves to the discussion of the results obtained in the case 
F = 1 A^^. A study of the half- filled case for fixed a and 
varying F has already been published ||3l|| . 

Given the operators Ai and Bj, defined at sites i and 
j, we evaluate the correlation function (CF) {AiBj) in 
the ground state \G) of the effective Hamiltonian H* ac- 
cording to: 

= (G|(Op^|e^e^A,S,e-^e-^|Op„)|G) . (12) 

The nature of the ground state is investigated for differ- 
ent fillings by evaluating the dependence on the lattice 
constant of the spin and charge CFs 



Xn 



IE 



(Si-Si+i) 



(13) 
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Xn 



(14) 
(15) 
(16) 



with Xn'' and x^n indicating the nearest and next-nearest 
neighbor CFs for the spin (s) and charge (c) channel, 
respectively. In addition, we have analyzed the behavior 
of the on-site, nearest- and next-nearest neighbor charge- 
deformation CFs, given respectively by: 



j 



(17) 
(18) 

(19) 



It is easy to show that the above charge-deformation CFs 
can be expressed as weighted sums of the charge CFs with 
the real-space displacement parameters acting as weights. 
In the case of dispersionless phonons their expressions are 



'J2a 



2hn 



2nn 



2hn 



[nSioc + 25aXn + 52aXnn) (20) 
[{Sloe + 52a)x'n + ^aXnu + ^Sa] (21) 
{25aXn + SlocXnn + "'52a) , (22) 



where n — N/L is the electron density and 5ioc^ 5a and 
are real-space displacement parameters, associated 
with the coupling of an electron to a lattice distortion, 
respectively on the same site, on nearest and on next- 
nearest neighbor sites. Recalling that for a four-site chain 
the wavevector q takes the 4 values —it /2a, 0, it /2a and 
7r/a, and assuming 5^^ /2a = 5^ /2a and q;_^/2q = a,r/2a, 
one gets that the expressions of 5ioc, 5a and 52a in terms 
of the g-dependent displacement parameters are 



Sloe = 5ii = ^ ((5o + 2(5jr/2a + 5^/a) 
5a = 5iA+l = - ((5o - (^Tr/a) 
52a = 5i^i+2 = J (^0 - 25jr/2a + ^7r/a) 



(23) 
(24) 
(25) 



(similar expressions are obtained for the real-space 
squeezing parameters ) . It is worth noting that in the 
case of a single electron (i.e. = 1), Eqs.(pO|)-(p^) give a 
simple relation of proportionality between the functions 
$ and the parameters 5 evaluated at the same site dis- 
tance. 

In the numerical analysis for N > 2 we have paid spe- 
cial attention to the effect on the phase diagram of the 



phonon-induced charge interaction, by distinguishing the 
results obtained by including both nearest {V*^^^) and 
next-nearest {V*^^2) inter-site charge interactions, from 
those obtained for ViiJ^2 = 0. The behavior of the corre- 
lation functions and the renormalized electronic interac- 
tions discussed in the next sections, always refers to the 
case of ViiJ^x and y*i±2 both non- vanishing, with the pa- 
rameter F chosen equal to 1 A~^. We have also verified 
that at any filling the parameters X*, P* and J* = J 
are very small and have no influence on the phase dia- 
gram. Therefore their behavior will not be shown in the 
following. We shall analyze the results for the fillings 
N ~ 1, 2, 3, 4 separately. 



A. The case iV = 1 

We treat this case mainly to test the reliability of our 
approach as compared with other methods developed in 
the literature ffl [7|, P,[T8| ,p6|. The phase diagram for 
= 1 is shown in Fig!^. 




FIG. 1. Phase diagram for iV = 1 and T = lA"^ The 
labels MP and TP stand for mobile polaron and trapped po- 
laron, respectively. 

The boundary separates two regions. In the lower one 
the polaron has a non- vanishing hopping amplitude, with 
the deformation involving, though very weakly, the first 
and second neighbor sites. In the upper region t* vanishes 
together with <^a and $2a, so that the polaron is trapped 
and the deformation is strictly local. This picture is in 



4 



substantial agreement with the results presented in Refs. 
~T1. 



B. N = 2, the quarter-filling case 

The phase diagram concerning the quarter-filling case 
{N = 2) is plotted in Fig]| for T = 1 A"! comparing the 
cases with and without V*^^2■ We can distinguish three 
main regions: a disordered (D) phase without magne- 
tization, a phase in which the two electrons behave as 
separated polarons (SP) in a superposition of configura- 
tions of the 1 — — 1 — type, and a phase characterized 
by the formation of an on-site bipolaron (OSB), i.e. a 
superposition of configurations of the 2 — — — type. 
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FIG. 2. Phase diagram for = 2 and P = 1 A"^ The re- 
gion labelled D, SP and OSB are characterized by the presence 
of, respectively, charge and spin disorder, separated polarons 
and on-site bipolaron. 

We see that in this case the Coulomb interaction be- 
tween next-nearest neighbor sites V*^J^2 affects the phase 
diagram rather strongly. Indeed, while there is no sig- 
nificant effect on the boundary of the OSB region, the 
inclusion of the effect of 1^* -1-2 leads to a considerable 
enlargement of the D region at the expenses of the SP 
region, additionally changing the D-SP transition from 
sharp to smooth. 
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FIG. 3. Spin (top panel), charge (middle panel) and 
charge-deformation (bottom panel) correlation functions for 
N = 2aX g'^/hQ. = 0.1 eV. 

Let us discuss in more detail the D-SP crossover for 
small polaron binding energy, choosing /Ml — 0.1 eV. 
As one can see from Fig.^ (top panel), the behavior of 
the spin CFs reveals a competition between different spin 
configurations in the ground state, without a clear pre- 
dominance of one of them. In particular Xnn almost 
zero while Xn is negative in the D region revealing a ten- 
dency towards the AF configuration. On the other hand, 
the charge CFs in the middle panel of Fig.|^ show that 
above a critical value of a (in the case considered in the 
figure equal approximately to 2.8 A), the charge distri- 
bution changes, giving rise to a pattern 1-0-1-0 
of alternating singly-occupied and empty sites which is 
characteristics of what we have called the SP region. 

On the contrary, in the D region, Xn and Xnn are 
smoothly changing with a, revealing the presence of a 



5 



fluctuating regime with disordered charge distribution. 
Let us stress that in the SP configuration there is no 
predominant character of either FM or AF correlations, 
since the shghtly negative values of Xnn cannot be con- 
sidered as a sign of the onset of an antiferromagnetic 
order. The displacement parameters and the charge- 
displacement correlation functions reveal additional in- 
teresting details of the transition. 

Actually, the bottom panel of Fig.|^ shows that there 
is a smooth, yet pronounced, change in the phononic 
regime for a approximately greater than 2 A, thus well 
inside the D region. Indeed, for small values of a, all the 
charge-deformation CFs are non-negligible and compara- 
ble. Then, when a is increased above 2 A , $/oc and $2a 
grow, both converging to the value 0.25, while <I>a de- 
creases becoming negligible above a = 3 A. The resulting 
smooth crossover from the D to the SP phase shows that 
already in the D phase there are precursory effects of the 
charge ordering of the 1 — 0—1 — type that character- 
izes the SP phase, stable at higher values of the intersite 
distance. 



pointing out that the D and SP states are very close in 
energy, so that even a small interaction like T^i,i±2 can 
significantly influence the position of the phase bound- 
ary. 

Considering values of the polaron binding energy 
/Ml greater than approximately 0.7 eV, we have that 
the increase of a induces a transition from the disordered 
phase to a new charge ordered phase, in which the local 
effective interaction U* becomes strongly attractive, thus 
leading to the formation of an on-site localized bipolaron 
(OSB). The nature of the OSB region and the driving 
mechanism leading to its formation are clarifled by the 
behaviour of the renormalized electronic parameters. As 
one can see from Fig.||, in the OSB phase the on-site 
Coulomb interaction is so strongly renormalized that it 
becomes attractive, while at the same time the kinetic 
energy vanishes. This behavior of U* and t* leads to the 
formation of a localized on-site bipolaron. 
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at g'^/hQ. = 0.9 eV 



FIG. 4. Renormalized electronic interaction parameters for 
A/' = 2 and g"^ /hO. = 0.1 eV. 

The corresponding behavior of the effective two-body 
parameters is shown in Fig.^. In particular, it is interest- 
ing to notice that a signature of the two phononic regimes 
in the D phase discussed above is clearly visible in the 
behavior of the next-nearest neighbor interaction V*^^2 j 
which already at a ~ 1.7 A turns from repulsive to at- 
tractive, while both U* and V*^J^■^ keep being repulsive. 
This situation clearly promotes the SP type of charge 
order, which, however, can be stable only when the effec- 
tive hopping amplitude t* becomes negligible. It is worth 



We have also found that in the OSB phase both Xn and 
Xnn vanish for any a, so that the phase is characterized 
by the absence of both spin correlations and local mo- 
ments. In addition, the charge correlation functions Xn 
and Xnn are also identically vanishing, while '^loc is the 
only charge-deformation CF which is different from zero. 
These results clearly indicate that in the OSB phase the 
charge, together with its accompanying deformation, is 
completely localized on one of the chain sites, in a su- 
perposition of configurations of the 2-0-0-0 type. 
This has also been verified by checking that in the ground 
state wavefunction all the components of the 1 — 1—0 — 
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type have negligibly small amplitudes. 

It is interesting to compare the present results with 
other studies on the behavior of two fermions coupled 
to lattice deformations. The main issue discussed in the 
literature is to clarify under which conditions a bound 
state of two polarons, i.e. a bipolaron, either on-site or 
inter-site, is stable, and possibly itinerant. In the large 
g'^/hfl region, dominated by the negative value of U* , our 
result of an essentially localized on-site bipolaron agrees 
with those in the literature We have also verified 

that, even imposing the vanishing of the bare inter-site 
charge repulsion Vij^ the OSB region is not changed. 

In the low g^/hfl region, we do not find mobile 
inter-site bipolarons, at variance with other studies 
|,U|Jl|,|l|,|l|,||. In assessing the stability of either 
type of bipolaron, it is very important to take into ac- 
count the bare contribution Vi^i±i to the effective inter- 
site charge interaction V*j . In the papers in which a stable 
inter-site bipolaron is found, this term, that is repulsive 
and usually rather large (a sizeable fraction of U), is how- 
ever neglected , or explicitly assumed to be 
renormalized to an overall attractive value Another 
recent paper emphasizes the physical relevance of the 
inter-site bipolarons. Within its LFA treatment of the 
long-range Frohlich electron-phonon interaction with dis- 
persive phonons, the bare Vij is included, but the equiva- 
lent of our V*j is estimated, on a purely phenomenological 
basis, to be very small, even between nearest neighbors. 

Our main conclusion on the formation of a bipolaron 
and, in particular, on its size, is that, when the complete 
set of bare two-particle interactions is taken into account, 
for N = 2 only on-site bipolarons are possible, this tak- 
ing place in the range of parameters corresponding to 
the OSB region of Fig.|. When stable, the OSB has a 
localized character, due to the simultaneous vanishing of 
the effective single-particle hopping t* and the effective 
pair hopping amplitude P*. Hence, our results do not 
support the presence of intersite mobile bipolarons for 
realistic values of the interactions, even if on the other 
hand they cannot definitely rule them out in infinite lat- 
tices. 

Besides the limitation of dealing with a cluster, one 
might also question the approach that we adopt here 
to evaluate Vij jl^, yielding rather high values of the 
bare intersite repulsion that, for physically reasonable 
values of the phonon frequency Q and of the coupling 
strength g, can be only partially compensated by the 
phonon renormalization. However, even if the absolute 
magnitudes of the various interactions we evaluate are 
model-dependent, their ratios should be of more general 
applicability. Indeed, it is known that in several classes 
of physical systems the nearest-neighbor Coulomb inter- 
action is a significative fraction of the on-site one. This 
is, for instance, the case of the high- Tc superconducting 
copper oxides, where theoretical estimates of Vi^i+i PJ33| 
give a repulsive value of the order of 0.2 eV. 



C. A'^ = 3, one-hole case 

The phase diagram for = 3 is shown in Fig.|[ We 
start by discussing the phases that are stable at low val- 
ues of the polaron binding energy, looking at the be- 
haviour of the correlation functions for g^ /Tift — 0.1 eV. 
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FIG. 6. Phase diagram for A = 3 and F = 1 A"^ The re- 
gions labelled NM, CO and FR are characterized, respectively, 
by vanishing magnetization, charge order and frustrated mag- 
netic order. 

As one can see from the upper panel of Fig.^, the be- 
havior of Xn shows little changes as a is varied and it is 
always negative. On the other hand, xtm positive and 
close to zero for a smaller than ~ 2 A and then it changes 
sign as a is increased, going to a saturation value of about 
—0.125. As far as the charge channel is concerned, we 
see from the middle panel of Figj^ that Xn exhibits slight 
variations with a, with a smooth increase concentrated 
in a range going approximately from 2 A to 3 A. On the 
other hand, up to a « 3A Xnn larger than Xni ^^^^ 
growing (up to a ~ 2 A) and then smoothly decreasing 
towards the saturation value 0.5 that both Xn ^'^d Xnn 
reach for a > 3A. 
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FIG. 7. Spin (top panel), charge (middle panel) and 
charge-deformation (bottom panel) correlation functions for 
N = 3 at g'^/nn = 0.1 eV. 



terizes the ground state of the system at higher values of 
g^/hft. Finally, the charge distribution tends to become 
uniform as the system, for a > 3 A, enters the region that 
we label as "frustrated magnetism" (FR). Indeed, there 
Xn and Xnn Converge to the same unsaturated negative 
value, which expresses a tendency of the spins to align 
antiferromagnetically both on nearest- and next-nearest 
neighbor sites. The FR zone is also characterized by a 
reduction in the extent of the deformation around each 
occupied site, as indicated by $2a becoming the smallest 
charge-deformation CF. 

The renormalized electronic parameters (not shown) 
vary smoothly with a. They never change sign, with 
the exception of the long-range charge interaction which, 
always very small, changes from positive to negative near 
the NM-FR boundary. 

Let us now consider larger values of the polaron bind- 
ing energy, investigating in particular that part of the 
phase diagram, corresponding to values oig^ /fiQ. between 
approximately 0.25 A and 0.7 A, where a new charge or- 
dered (CO) phase develops between the NM and the FR 
regions for intermediate values of a. This will be done 
with a special attention to the character of the CO phase, 
by analyzing the behavior of the CFs at the representa- 
tive value /fi^l = 0.4 eV. As one can see from the up- 
per panel of FigJ^, for values of a falling in a range going 
approximately from 2.2 A to 3.7 A the spin CF Xn a'^^ 
Xnn S''^^ zero, indicating the absence of local moments 
and magnetic correlations in the ground state. On the 
other hand, by looking at the behaviour of Xn Xnn 
in the same range of values of a, we can see that charge 
correlations are vanishing on nearest-neighbor sites and 
tend to the maximum value compatible with the filling 
on next-nearest neighbor sites. This is a clear indication 
of the development of an ordered configuration where the 
charge is fully distributed on next-nearest neighbor sites. 



As far as the charge-deformation CFs are concerned, 
we see from the bottom panel of Fig.|^ that $;oc is the 
largest at any a, with <i>2a manifesting a significative in- 
crease exactly in the range of intermediate values of a 
in which Xnn reach its maximum. In addition $a turns 
out to be the smallest charge-deformation CF at any a, 
becoming equal to <I>2a only above a « 3 A. 

Combining the results for the CFs we conclude that 
for a lower than « 2A in the non-magnetic (NM) re- 
gion, the charge is distributed with fluctuating config- 
urations of almost empty and almost doubly occupied 
next-nearest neighbor sites. This causes both Xnn ^-ud 
Xnn to be negligible. Then, as a is increased beyond 2 A 
up to a w 3 A, there is, consistently with the behavior of 
^■20, a crossover region which can be seen as a precursor 
of the charge order of the 2 — — 1 — type, which charac- 
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FIG. 8. Spin (top panel), charge (middle panel) and 
charge-deformation (bottom panel) correlation functions for 
N = 3 at g'^/nn = 0.4 eV. 

The same general trend can be inferred from the 
charge-deformation CFs, shown in the lower panel of 
Fig.^. We see that the correlations between the charge 
and the deformations in the CO region become the 
strongest between second neighbor sites, with the first 
neighbor one being the smallest. This can be interpreted 
as a further indication that the ground state is given by 
a superposition of the charge configurations 1 — — 2 — 
and 2 — — 1 — 0, together with the symmetrically related 
0-1-0-2 and 0-2-0-1. This degeneracy of the 
ground state, related to the fact that we cannot have a re- 
ally broken-symmetry state in finite-size systems, is also 
the reason why in the CO region Xmi equal to 1 , rather 
than being equal to 2. The whole picture is also con- 
sistent with the behavior of the renormalized electronic 
parameters, featuring a non-negligible negative value of 



V^* -1-2 both in the CO and in the PR region. 

The filling iV = 3 can be viewed as the case of one 
hole in an half-filled configuration. There are two differ- 
ent predictions on the related ground state in the limit 
of large U* /t*. According to Nagaoka's theorem, a ferro- 
magnetic (FM) state should set in when [/ — > cxd with t 
remaining finite, while according to Refs. [^,||jl^,|l^ the 
hole would be a polaron in an antiferromagnetic (AF) 
background, where its phonon-induced localization would 
take place at a lower value of /TiVL due to concurrent 
effect of the ordered spin arrangement. In our approach 
we find that for iV = 3 none of the two magnetic ground 
states mentioned above is stable. On the contrary, in 
what we have called the frustrated region no magnetic 
effect prevails, this region being characterized by FM and 
AF correlations of comparable amplitude. While one can 
only speculate on the fact that we do not find an AF 
ground state, guessing that it is a finite-size effect, on 
the other hand the reason why we find no FM order is 
that the conditions under which Nagaoka's theorem holds 
are never satisfied. Indeed there is no parameter regime 
where the phononic renormalization yields U* 3> , with 
t* remaining at the same time finite. 

D. AT = 4, the half-filled case 

In Fig.0 we report the phase diagram at half-filling for 
r — lA , both with and without the phonon-induced 
next-nearest neighbor interaction. For intermediate to 
large a we see that an antiferromagnetic (AF) phase is 
stable for small /hVl while a charge ordered state (CO) 
of the 2-0-2-0 type is stable for large jhKl. For 
small a and irrespective of g'^ /Kl, a non-magnetic (NM) 
ground state characterized by the absence of any spin 
correlations is found to be stable, wedging in between 
the AF and CO regions at small g^/hil, up to a ~ 2.5 A. 

We can also see that taking into account the effect of 
Vj* -1-2 leads to quantitatively significant modifications of 
the phase diagram. Indeed, even if the negative value of 
^ii±2 is rather small for low g^ /Ml values, it causes an 
appreciable expansion of the NM region towards higher 
g'^/hQ and a values, with a corresponding reduction of 
the CO and AF regions, due to the fact that for increasing 
g'^/hil (or a) V*^j^2 becomes negative, and thus opposite 
in sign to V*^j^^. 
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1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 
a (A) 

FIG. 9. Phase diagram for = 4 and T = 1 The re- 

gions labeUed NM, CO and AF are characterized, respectively, 
by vanishing magnetization, charge order and antiferromag- 
netism. 

Let US now discuss how the system responds as the 
value of a is increased, for /Ml = 0.1 eV. As shown in 
Fig.^, the system passes from the NM to the AF regime. 
The position and the nature of this transition can be eas- 
ily identified by looking at the behavior of the correlation 
functions plotted in Fig.O. 
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FIG. 10. Spin (top panel), charge (middle panel) and 
charge-deformation (bottom panel) correlation functions for 
iV = 4 at g'^/hQ. = 0.1 eV. 

For a lower than ~ 2.2 A, Xn Xnn ^'"^ close to zero 
and almost coinciding, which corresponds to a regime 
without significative magnetic correlations. On the other 
hand, for higher values of a one has Xn ~ —0.5, which cor- 
responds to an antiparallel configuration of spins on first 
neighbors, and Xnn ~ 1/4 which indicates a configura- 
tion of parallel spins on next-nearest neighbors. This be- 
haviour of Xn and xtm clearly shows that the ground state 
has antiferromagnetic correlations. The middle panel of 
the Figjl^ gives indications on the charge distribution. In 
the NM region we can see that even though the site occu- 
pation is still fiuctuating because of the itinerancy of the 
carriers, nonetheless there is some precursor of the CO 
phase, as a consequence of the fact that the values of Xnn 
are considerably higher than those of Xn- Conversely, in 
the AF region Xn and Xnn become constant and equal to 
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one, thus indicating that there is one electron per site on 
the average. 

The same quahtative trend is found in the charge- 
deformation CFs, shown in the bottom panel of Fig.[Tc|. 
At small a the charge-related deformation appreciably 
affects also first and second neighbor sites. Then, as 
the NM-AF boundary is approached, ^loc tends to stay 
constant, while $20 grows considerably and <i>a becomes 
small, resulting in a behaviour that can be considered as 
a precursor of the CO state. When the AF state sets in, 
all three charge- deformation CFs come to coincide at a 
rather large value, indicating that, in a regular pattern 
of site occupancy, all orders of neighbors are equally af- 
fected. In terms of the size of the polarons, one can say 
that for iV = 4 they tend to show an extended character, 
i.e. the associated deformation develops beyond the site 
where the charge is located. 

The behavior of the renormalized electronic parame- 
ters (not shown) is such that in the AF region t* is van- 
ishingly small, driving the system in a regime of electron 
localization where the mechanism of superexchange dom- 
inates. Conversely, in the NM region, the hopping term is 
only weakly renormalized, because the localization effect 
due to the displacement is counteracted by the squeezing, 
so that the electrons are rather mobile. Moreover, U* and 
V*^^i are comparable in magnitude, thus preventing the 
formation of local pairs as well as of long range ordering. 
In terms of an equivalent t — J model, the NM ground 
state does not show any well-defined spin or charge or- 
der since t* remains large enough to have two concurrent 
effects: on one side it opposes a complete charge dis- 
proportionation, and on the other side it gives rise to a 
kinetic exchange proportional to 2t*^ /{U* — V*^i) strong 
enough to allow for some magnetic coupling effects. As in 
the previous cases, the energy difference between the NM 
and the AF states is small enough for the tiny negative 
value of V^i* ±2 to have a sizeable effect on the position of 
the phase boundaries. 
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FIG. 11. Spin (top panel), charge (middle panel) and 
charge-deformation (bottom panel) correlation functions for 
A/' = 4 at g^/nn = 0.4 eV. 



When g^/fin — 0.4 eV, a reentrant behavior can be de- 
duced from the phase diagram in Fig.^ as a is increased. 
This is made evident by the behavior of the correla- 
tion functions in Fig.^. For a between approximately 
2 A and 3.6 A, 

Xn is equal to zero (middle panel), while 
Xnn assumes the largest value (~ 4), clearly indicating 
a CO configuration of alternating empty and doubly oc- 
cupied sites. As expected, in this same range of a the 
spin correlation functions (upper panel) are zero since 
configurations with zero and two electrons on each site 
prevent the formation of magnetic moments. As far as 
the charge-deformation CFs are concerned, we see that 
in the CDW region the largest values are taken by <I>2a, 
with the nearest-neighbor CF remaining always very 
small, and the on-site CF ^loc taking values which are 
intermediate between $20 and ^loc- This can clearly be 
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considered as a further evidence of the formation of a 
CO state of the 2 — — 2 — type. When a Ues out of 
this intermediate range, the correlation functions behave 
essentially as in the case /Ml = 0.1 eV. 

The renormalized electronic parameters (not shown) 
fully support this picture. In the CDW region the elec- 
tron mobility is equal to zero, V*^^2 is attractive and 
is repulsive and comparable in magnitude to U* . 
This evidently tends to promote the configurations with 
alternating empty and doubly occupied sites. To the ex- 
tent that the effects of including V*j can be neglected, 
our results qualitatively agree with other studies |^ . 

We end this Section by noting that the discontinuities 
in the correlation functions and the renormalized param- 
eters seen in the above figures, should not be interpreted 
as a signature of a real first-order localization transition. 
They are rather a consequence of a drawback affecting 
all variational methods, which, as it is well-known, can 
only provide an upper bound to the true ground state 
energy. This leads to a violation of the exact results 
demonstrated by Lowen Q , according to which the self- 
trapping transition generated by an interacting term of 
the Holstein type is an analytical crossover for any value 
of the electron-phonon coupling, and not an abrupt (non- 
analytical) phase transition. 

IV. CONCLUSIONS 

The method that we have used to obtain an effec- 
tive polaronic Hamiltonian from the generalized Holstein- 
Hubbard Hamiltonian, based on the use of the displace- 
ment and squeezing transformations, while allowing tech- 
nically simple calculations, yet permits to grasp some 
basic physical effects, such as, for instance, the phonon- 
induced renormalization of the electronic interaction pa- 
rameters. In this respect, this approach might be of in- 
terest, for instance, in the theoretical analysis of high- 
Tc superconductivity, where it can provide clues to se- 
lect among electron-phonon models or purely electronic 
models. If, in the comparison between theory and exper- 
iments, one refers to phonon-rcnormalized, and not bare, 
interactions, then even phenomena such as the isotope 
effect could find an explanation within the framework 
of a model with only electronic degrees of freedom. An 
improvement of the approach followed here may be ob- 
tained by suitably dealing with two limitations, that is, 
(i) the factorization of the fcrmion-boson wave function 
into separate fermionic and bosonic factors, and (ii) the 
use of a squeezed phonon wavcfunction (whose intrin- 
sic limitations have been discussed in Ref. [0) for the 
elimination of the phononic degrees of freedom. Further 
studies in this direction are planned for the next future. 

The results we have presented support our claim that, 
in working out an effective Hamiltonian through displace- 
ment and squeezing transformations, it is important to 



take into account the momentum dependence of the dis- 
placement and squeezing parameters, to avoid losing a 
significant part of the physical content of the model. 
We have indeed shown that the effective phonon-induced 
inter-site charge interaction that one gets [ p5p^ ] by go- 
ing beyond LFA, should be carefully taken into account, 
as it has appreciable effects on the phase diagram, par- 
ticularly for = 4 and N — 2. This latter case is par- 
ticularly significant, because theoretical studies of 
the occurrence of doping-induccd inhomogeneities in the 
charge distribution in both colossal magnetoresistance 
manganites and superconducting cuprates, point to the 
importance of the long-range charge interaction, but have 
difficulties in justifying its presence in metallic systems 
where it should vanish beyond a very short distance, due 
to the screening induced by the itinerant carriers. Dif- 
ferently from the purely electronic Vij, the effective V*j 
can be seen as an independent parameter providing the 
needed long-range interaction through a mechanism de- 
pending basically on the phonons, and only weakly on 
the doping. In real systems such as the superconducting 
cuprates this mechanism can be ascribed to the Holstein- 
type electron-phonon coupling generated by the ion dy- 
namics out of the conducting layers. Additionally, V*^^2 
can vary both in amplitude and, most importantly, in 
sign, depending on the phononic state and on the value 
of the lattice constant. 

lYvoTn a more general point of view, our results show 
that there exist, for Af = 2, 3, regions of the phase dia- 
gram where the ground state of a hole-doped correlated 
system with strong electron-phonon interaction, on one 
side has appreciable coexisting spin and charge correla- 
tions, while on the other side is such that the carriers 
retain, to a large extent, their itinerancy. This suggests 
that the debate on whether charge or spin fluctuations 
in the normal state are responsible for high-temperature 
superconductivity might be viewed from a different per- 
spective, namely, by considering that, under appropriate 
conditions, both fluctuations can coexist, possibly coop- 
erating to build the superconductive phase. 

A significant result is that, for A^ > 2, in the dis- 
ordered regions of the phase diagrams for low values 
of a and /TiVL, we find wide subregions characterized 
by strong precursory charge-ordering effects, whose am- 
plitude can change very rapidly with the lattice pa- 
rameter. The effect is particularly evident in the case 
N = 2. This might be of interest with respect to the 
recently reported |2|] microscopically coexisting charge 
and structural dishomogeneities in underdoped and op- 
timally doped La2-a;Srj;Cu04 samples. 

Finally, we have been able to model the effects of a 
varying lattice spacing. This might be of relevance in 
all the cases in which compression (or expansion) of the 
lattice causes sharp phase transitions. In particular, our 
results show that in the range of a where a given state 
is stable, even if that state does not show a fully devel- 
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oped symmetry-breaking order parameter (CO or AF), 
still the behavior of some correlation functions may be 
characterized by drastic changes as the site distance is 
varied. 
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